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We compute the structure of a quantized vortex line in 
a harmonically trapped dilute atomic Bose-Einstein conden- 
sate using the Popov version of the Hartree-Fock-Bogoliubov 
mean-field theory. The vortex is shown to be (meta)stable in a 
nonrotating trap even in the zero-temperature limit, thus con- 
firming that weak particle interactions induce the condensed 
gas a fundamental property characterizing "classical" super- 
fluids. We present the structure of the vortex at ultralow 
temperatures and discuss the crucial effect of the thermal gas 
component to its energetic stability. 



PACS number(s): 03.75.Fi, 05.30. Jp, 32.80. Pj, 67.40. Db 

Dilute atomic Bose-Einstein condensates (BECs) have 
been the subject of numerous experimental and theoret- 
ical studies since the first landmark experiments in 1995 
These quantum fluids are unique in providing an 
opportunity to investigate how the phenomenon of Bose- 
Einstein condensation is affected by weak particle inter- 
actions. The fundamental question, whether superflu- 
idity can be sustained by particle interactions in such 
systems, has still remained partly open despite extensive 
research. 

The existence of stable, dissipationless vortices is a 
feature characteristic to superfluid behavior. This as- 
pect on the superfluid properties of dilute BECs is espe- 
cially timely due to the recent experimental realizations 
of vortices in such systems [^,||. The stability of vortex 
structures has been under vigorous theoretical analysis 
P'p^. It has been shown that vortices in cylindrically 
trapped condensates are unstable within the Bogoliubov 
mean-field approximation P- p^Jl^ , p^ unless the system 
is continuously driven by a suitable rotating perturbation 
pO| , p^p^ . On the other hand, by taking into account ef- 
fects of the thermal gas component in the system, vortices 
have been shown to become (meta) stable at high enough 
temperatures even in a nonrotating trap, or when a suit- 
able external pinning potential is applied Q . 

The absence of dissipation in the superfluid flow im- 
plies that a circulating current persists even when the 
system is not rotated by an external perturbation. "Clas- 
sical" superfluid behavior thus implies stability of vor- 
tices even in a nonrotating vessel. Strictly speaking, 
such states are local minima of free energy and only 
metastable, the global minimum naturally corresponding 



to a nonrotating state. However, in order to investigate 
the superfluid properties of dilute boson condensates, it 
is just this kind of local energetic stability of vortices that 
needs to be clarified. 

In this Letter we present results of computations con- 
cerning the structure of a cylindrically trapped dilute 
BEC containing a vortex line, and show that such states 
are locally energetically stable at all temperatures T < Tc 
(Tc denotes the critical temperature of condensation) 
even in a nonrotating trap without additional pinning 
potentials. In this respect such systems are indeed shown 
to behave like "classical" superfluids. We also discuss the 
obvious discrepancy between this result and the predic- 
tions of the Bogoliubov approximation, and the role of 
the thermal gas fraction in stabilizing the vortex state in 
the zero-temperature limit. 

We consider a dilute Bose-condensed gas consisting 
of atoms with mass m, trapped by a radial harmonic 
potential T4rap(r) — ^mui^r^ in cylindrical coordinates 
r — (r, 6*, z). The particle interactions are modelled by an 
effective low-temperature contact potential T^nt(r,r') = 
g5{v — r'), with the bare interaction constant g and s- 
wave scattering length for binary collisions a related by 
g — Ai:h^a/m. Assuming the total particle number in the 
condensate to be large enough to justify grand canonical 
formalism, the equilibrium condensate wavefunction 0(r) 
satisfies the generalized Gross-Pitaevskii (GP) equation 



{Ho + t/c(r)|^(r)|2 + 2(7e(r)p(r)}0(r) = ^lcj>(v), (1) 

where Ho — — ft^V^/2m -I- Vtiap(r) is the bare single- 
particle Hamiltonian for the trap, the chemical poten- 
tial, and p(r) the density of the noncondensed gas. Func- 
tions C/c(r), Uc{y) are effective interaction couplings for 
various mean-field approximations — for the Popov ver- 
sion they are simply chosen as Uc{y) = Ue{r) = g. In ad- 
dition, we have performed computations within so-called 
Gl and G2 approximations, which are gapless mean- 
field theories taking into account effects of the back- 
ground gas on colliding atoms, neglected in the Popov 
approximation j2^j2^. In the Gl version Uc{r) = g and 
Uc{r) — g[l + A(r)/(/)^(r)], where A(r) is the anomalous 
average of two Bose field operators; for G2 one chooses 
Uc{r) ^ Uc{r) = g[l + A{r)/(l)'^{r)]. The usual procedure 
of diagonalizing the mean-field Hamiltonian by a Bogoli- 
ubov transformation yields coupled eigenvalue equations 
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for the bosonic quasiparticle amplitudes Uq{r), Vg{r), and 
the eigenenergies Eg of the form 

Cug{r) + C/e(r)02(r)v^(r) = EgUg{r), (2a) 
Cvg{r) + [/c(r)0*'(r)M,(r) = -EgVgir). (2b) 

Above, C = Ho - H + 2Uc{r)\(f>{r)\^ + 2t/c(r)p(r) and 
q denotes quantum numbers specifying the quasiparticle 
states. In addition, we have selfconsistency relations for 
the noncondensate density p(r) and for the anomalous 
average A(r): 

p(r) = ^[(|^,(r)p + \vgir)\')niEg) + \vgir)\% (3) 
A(r) = Y.i^Ugir)v;ir)niEg) + Ug{r)v;{r)], (4) 

where n{Eq) = {e^t/^^'^ — is the Bose distribution 
function. The anomalous average is ultraviolet divergent, 
and we renormalize it by subtracting the last term in the 
sum of Eq. (|) |2|]. 

Considering a condensate penetrated by a single vortex 
line, we search for solutions of the form 0(r) = 0(r)e*™^, 
where m denotes the number of circulation quanta of the 
vortex. In this Letter we restrict to the case m = 1 due to 
instability of multiquantum vortices [0,^6|jl^. By uti- 
lizing the cylindrical symmetry, Eqs. (|ip and (||) can be 
reduced to radial equations, which we discretize using a 
finite-difference method. Dirichlet boundary conditions 
are imposed at r = i?, the radius R chosen large enough 
for finite size effects for the structure of the vortex to be 
negligible. In the z-direction we impose periodic bound- 
ary conditions at z = ±i/2, thus modeling a system in 
the limit of a very weak axial trapping potential. Due to 
cylindrical symmetry, the quasiparticle amplitudes can 
be chosen to be of the form 



Uq{v) = y^(^)e»9.(2Vi)-+»(?«+™)9^ 



(5a) 
(5b) 



where qe and qz are integer angular and axial momentum 
quantum numbers, respectively. Discretization trans- 
forms Eqs. (H) to a narrow-banded matrix eigenvalue 
problem, which we solve using the Lanczos algorithm im- 
plemented in the ARPACK subroutine libraries p^ , p5| . 
The nonlinear Gross-Pitaevskii equation is solved us- 
ing finite-difference discretization and an overrelaxation 
method. For a given value of the chemical potential /i, 
the solution of the GP equation and the noncondensate 
density are integrated to find out the total particle num- 
ber, and the process is iterated until the chemical poten- 
tial corresponds to the preset total number of particles. 
The solution of the GP equation can be mapped to a 
zero-energy solution of Eqs. (|^), thus providing a test for 
the accuracy and consistency of the numerical methods 
used to solve these equations. We search selfconsistent 



solutions for Eqs. (|l|)-(^ using an iterative scheme: The 
condensate wavefunction and chemical potential corre- 
sponding to a preassigned total number of particles are 
computed by solving the GP equation. Using the quasi- 
particle states obtained from Eqs. (^, new mean-field 
potentials are computed using the selfconsistency equa- 
tions. The whole procedure is repeated until convergence 
to a desired accuracy. In order to stabilize the iteration 
at low temperatures, we use underrelaxation in updating 
the mean-field potentials. 

The physical parameter values for the gas and the trap 
were chosen to be the same as in Ref. 13 1. The gas 
consists of sodium atoms with mass m = 3.81 x 10~^^ kg 
and s-wave scattering length a = 2.75 nm. The frequency 
of the trap is Vr — LOr/lir — 200 Hz, and the density of the 
system is determined by treating N = 2 x 10^ atoms in 
the computational domain with dimensions i? = 20 /xm 
and i = 10 /xm. The critical condensation temperature 
for the system is approximately Tc f» 1 /xK. The spatial 
grid for discretizing the quasiparticle eigenequations was 
chosen dense enough to guarantee a relative accuracy of 
10"'^ for the energy of the lowest excitation, and better 
than 10"'* for the other states. 

Results of the selfconsistent computations performed 
within the Popov approximation are shown in Figs. 1-4. 
Figure |l| displays part of the selfconsistent quasiparticle 
excitation spectrum for the condensate vortex state at 
the temperature T — 100 nK. Only states with qz — 0, 1, 
which contain the excitations of lowest energy, are shown. 
These states determine the local energetic stability of the 
vortex configuration: If there exists a quasiparticle exci- 
tation with negative energy, the condensate can lower its 
energy by exciting this negative bosonic mode, the state 
thus being prone to collapse. For single-quantum vor- 
tices the lowest excitations with qg = —1, the so-called 
Kelvin mode states, form standing waves localized at 
the vortex core, the energies of these states being sub- 
stantially lower than those for other angular momentum 
quantum numbers. The sign of the energy of the lowest 
Kelvin mode excitation, the lowest core localized state 
(LCLS), determines the local energetic stability of the 
vortex. The LCLS is denoted by a diamond in Fig. 
Due to the very low energy of the LCLS compared to 
the other excitations, its contribution to the mean-field 
potentials in the core region becomes dominant in the 
zero-temperature limit. On the other hand, even minor 
changes in the potential functions affect the energy of 
the LCLS, thus altering its contribution drastically via 
the sensitive Bose factor. This behavior makes the itera- 
tion process extremely delicate at low temperatures [ p^ , 
requiring the use of suitable underrelaxation methods. 

The temperature dependence of the excitation energy 
of the LCLS is shown in Fig. |[ By carefully adjusting 
underrelaxation in the iteration procedure, we are able to 
compute the selfconsistent structure of the vortex down 
to the temperature of T = 1 nK. The convergence at 
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these temperatures is slow, but unquestioned. Especially, 
it was confirmed that the iteration converges to the same 
solution irrespective of the initial Ansatz used. However, 
it seems obvious that with additional computational ef- 
fort one could reach even lower temperatures. The result 
shows that although the energy of the LCLS approaches 
zero in the low-temperature limit, as expected, the vortex 
configuration is locally stable at all temperatures T < Tc 
I p^ . We obtained the same result also using the Gl and 
G2 approximations. 

The local energetic stability of the vortex within the 
Popov, Gl and G2 ||2^ approximations is to be com- 
pared with the instability predicted by the Bogoliubov 
approximation |^ , which does not treat the thermal part 
of the gas selfconsistently. Figure || shows the density 
profiles of the condensate and the noncondensate at ul- 
tralow temperatures. The thermal gas is concentrated in 
the vortex core, where it fills the space left by the con- 
densate and exerts an outward pressure on it, preventing 
the condensate from collapsing into the core. Selfconsis- 
tent treatment of the thermal gas fraction is thus crucial 
in determining the stability of vortices in such systems, 
as seen in the discrepancy between the predictions of the 
Popov, Gl and G2, and, on the other hand, the Bogoli- 
ubov approximation. At temperatures of the order of 
Tc, the validity of the Bogoliubov approximation is ex- 
pected to be questionable due to the substantial thermal 
gas fraction. It is to be noted that it can fail also in the 
low-temperature limit in certain respects, the stability of 
the vortex state being one example. This is due to the 
nonvanishing, residual noncondensate fraction present in 
interacting systems even in the zero-temperature limit. 
Figure ^ displays the computed thermal gas fraction at 
ultralow temperatures. It clearly shows the residual ther- 
mal fraction, which stabilizes the vortex at ultralow tem- 
peratures. 

It is to be noted that the positivity of the LCLS energy 
implies that vortices slightly displaced from the symme- 
try axis of the trap precess in the direction opposite to 
the condensate flow around the core j2^]. Recent exper- 
iments, however, show precession in the direction of the 
flow, except for a minority of the so-called rogue vor- 
tices [Q. This seems to imply a negative LCLS energy, 
in agreement with the zero-temperature Bogoliubov ap- 
proximation. The apparent discrepancy between the ex- 
periments and the predictions of the selfconsistent ap- 
proximations could be due to insufficient thermalization 
in the experiments of the gas in the region of the (mov- 
ing) vortex core. Wc expect that investigation of the 
validity of this assumption will further clarify the role of 
quasiparticles in stabilizing vortices in weakly interacting 
Bose-Einstein condensates. 

In conclusion, we have computed within selfconsis- 
tent mean-field theories the structure of a cylindrically 
trapped, dilute atomic Bose-Einstein condensate pene- 
trated by a vortex line. The vortex state is shown to be 



locally energetically stable in a nonrotating system even 
in the zero-temperature limit, thus confirming the system 
to act in this respect like a superfluid. The thermal gas 
concentrated in the vortex core is shown to have a cru- 
cial effect in stabilizing the vortex state even at ultralow 
temperatures, due to a residual noncondensate fraction 
present in interacting systems in the zero-temperature 
limit. 
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FIG. 1. Part of the selfconsistent quasiparticle excitation 
spectrum at temperature T = 100 uK for the harmonically 
trapped condensate containing a single-quantum vortex line. 
Only states with axial momentum quantum numbers qz = G 
(dots) and qz = 1 (circles) are shown. The lowest Kelvin 
mode state, the LCLS, is denoted by a diamond. The inset 
presents a blow-up of the lowest part of the full spectrum. 

FIG. 2. Energy of the lowest excitation, the LCLS, at tem- 
peratures T = 1-200 nK (Tc « 1 ^iK). The finite positive 
value of the lowest excitation frequency implies the vortex to 
be locally energetically stable. 

FIG. 3. Condensate and thermal gas density profiles at 
temperatures T = 1 nK, 100 uK, and 200 uK. The noncon- 
densate fills the space left by the condensate at the vortex 
core, thus exerting an outward pressure which stabilizes the 
vortex configuration. 



FIG. 4. Noncondensate fraction as a function of temper- 
ature for the vortex configuration. Due to particle interac- 
tions, there exists a residual thermal gas fraction even in the 
zero-temperature limit. 
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